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Abstract 

This paper covers the review, status, and projected future of direct numerical simulation 
(DNS) methodology relative to the state-of-the-art in computer technology, numerical 
methods, and the trends in fundamental research programs. 

1. INTRODUCTION 

Most review articles tend to present an exhaustive list of research publications doc- 
umenting the status of a technology. Review papers by Kleiser and Zang 1 and Reed 2 list 
such a treasure of literature on direct numerical simulation (DNS) applied to boundary- 
layer transition. In addition, discussions relating to issues for compressible DNS 3-6 and 
the use of DNS for turbulent flows 7,8 are available to the interested reader. In this short 
essay, it would be more valuable to look at DNS in the light of required technologies than 
to rehash or condense the available literature within the present page limitations. Hence, 
this paper will take a different twist on DNS of the past and look at future computations 
subjective to the author’s views. 

2. COMPUTATIONAL FLUID DYNAMICS 

Computational Fluid Dynamics, usually abbreviated as CFD, has been used to refer 
to the computational solution of fluid dynamics problems. In 1976, Chapman 9 showed 
a progression of CFD capability with time (Table 1). Using computer-speed forecasts, 
viscous time-averaged and time-dependent CFD were projected to be possible during the 
late-1970s and mid-1980s. Here, we are in the late-1990s and both the viscous time- 
averaged solutions are possible but some controversy surrounds accuracy and interpreta- 
tion of the results. Time- dependent solutions are now possible on simplified geometries 
and at low Reynolds numbers. Also because CFD is an all-encompassing term (or “buzz- 
word”) which has been used to cover fluid mechanics from linear ordinary differential 
equations through full unsteady nonlinear Navier-Stokes equations, a hint of skepticism 
naturally arises when discussing the value of CFD results. This skepticism (at times) has 
led to the humorous Color Fluid Dynamics and has implied some mystical skills are re- 
quired to obtain CFD results. In some sense, significant experience is required to interpret 
the CFD obtained solutions, lending to the assumptions/limitations of any given code. To 
help clarify the adequacy of a solution, perhaps every future paper or presentation which 
makes use of CFD (or a code name, such as CFL3D, TLNS3D, etc) should not only state 
that the Navier-Stokes equations are not solved directly but rather they are modeled and 
highlight the limitations of the results due to the assumptions within the methodology 
(e.g., low-order methodologies utilized, potential turbulence model limitations, inadequate 
grid resolution, etc.). Hence, we should change the paradigm of our engineers and man- 
agers from CFD implying a definitive (or exact) quantitative solution of the fluid and 
aerodynamic problem toward CFD as a tool which models the aerodynamics, which is 
ever striving for quantitative accuracy. 

3. DIRECT NUMERICAL SIMULATION 

Direct Numerical Simulation or Direct Navier Stokes (DNS) is a CFD technology 
implying a nearly exact solution to a unsteady, nonlinear governing systems of equations. 
In aerodynamics, DNS is associated with a large-scale computationally intensive solution 
procedure which may consume hundreds to thousands of Cray Super-computing resources. 



The earliest use of DNS began in the 1970’s. Table 2 shows the progression of DNS usage 
through 1997. A computer search of the literature indicates that an exponential growth 
in the use of DNS as a CFD methodology has occurred. Like conventional CFD, a state 
of overuse of terminology may be appearing in the literature. Although the scientific 
community can gleam the limitations inherent in the use of DNS in a given study, de- 
cision making managers or those reviewing manuscripts may not be familiar with such 
limitations. Furthermore, because DNS has become equivalent to “exact solution,” dis- 
crepancies in solutions, or comparisons, may be incorrectly attributed to behavior of the 
aerodynamics (flow physics) or to faults in experiments or theory. 

Table 1. Status of computational aerodynamics. 9 
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Table 2. Use of Direct Numerical Simulation (DNS) in the literature. 


Year 

Citations 

before 1970 

0 

1970-79 

12 

1980-89 

116 

1990-97 

536 


Perhaps the CFD and DNS community have suffered because the terminology has 
lacked sufficient definition to discriminate between the various levels of approximation of 
the governing equations. In some cases, DNS has been used for studies where the temporal 
and spatial resolutions have been insufficient to capture the quantitative behavior of the 
fluid dynamics. Hence, the sub-grid scale structures have been modeled in an ad hoc 
manner. Unresolved solutions often have value to study qualitative features of the problem 
in a less costly manner compared to a fully resolved study. Perhaps, we need to agree 
on what DNS implies and when DNS becomes an ad hoc LES procedure, and what we 
term conventional CFD. Perhaps additional definitions or conventions for the use of these 
terms and the introduction of some new abbreviations are in order. Some of these may 
take the form Course-Grid DNS (CGDNS) or ad-hoc LES (AHLES) or steady, low-order 
CFD (SLOCFD). 

Returning to the subject matter — past, present and future of DNS, the technical 
challenges associated with DNS involve extreme memory and computational speed re- 





quirements. Concerning memory, 3D simulations of solar convection, for example, could 
require on the order of 10 30 grid points. 10 Simulations of a planetary boundary layer 
could require 10 18 grid points, 10 while computers of the 1990s can handle only about 10° 
grid points. Concerning performance, DNS could be attainable for flow past an airfoil 
(Re ~ 10 8 ) if teraflop (10 12 flops) performance wave available. With exaflop (10 18 flops) 
performance, DNS could be enabled for flow past a complete aircraft 11 . Extending the use 
of DNS to more complex configurations is significantly linked to computer advancements 
and/or novel faster, and more accurate algorithms. And clearly, because LES is a fraction 
of the cost of DNS, large-eddy simulations have the potential to make large contribu- 
tions to the aerospace industry before full scale DNS is possible. Hence, my remaining 
comments apply to LES as well as DNS. 

4. COMPUTER ADVANCES 

From the 1976 presentation by Chapman, 9 computer speed versus year the computer 
was available is reproduced in our figure 1 and amended to include additional Cray prod- 
ucts and a projected computer with pentaflop performance in 2010. In relation to the 
predicted performance gains of 81%/year and 112%/year, computer technology advances 
have slowed. 
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Figure 1. Computer speed versus year available. 9 

In addition to raw computer speed advances, numerous multi-processor, parallel com- 
puters have entered the market for potential use in solving engineering and science prob- 
lems. Intel IPSC/860 and Paragon, IBM SP1 and SP2, and various other computers have 
been used to examine the feasibility of CFD on parallel architectures. Shown in figure 
2 are the cost of major kernels of a spatial DNS code, 12 the communication cost, and 
the total computing cost on the IBM SP2. Clearly, the SP2 becomes competitive with 
the same code using one processor of a Cray C90, in spite of the communication penalty 
incurred with a parallel computer. 

Also shown in figure 2, a single 3D computational domain can be conveniently (and 
efficiently) split among the processors by dividing the domain equally among the comput- 
ing nodes (division among 4 nodes is shown for example in figure 2). So parallel computers 



potentially buy performance gains and additional memory for the application of DNS to 
larger problems. For example, using a parallel computer, teraflop (10 12 flops) performance 
could be attainable and should enable DNS of flow past an airfoil. 




Figure 2. Domain decomposition for parallel computing and computer cost versus number 
of processors. 12 

5. NUMERICAL/ ALGORITHMIC METHODS 

While the speed increases bought by faster computers goes significantly to what can 
be accomplished using DNS, smarter programming, faster algorithms, and novel theoret- 
ical tools should receive continued emphasis to make progress in DNS. 

For periodically assumed flows, fast Fourier series methods have enabled numer- 
ous temporal DNS studies. For the spatial DNS approach, high order (>4th) finite- 
difference methods are commonly used in DNS codes. Due to the advances made by 
Lele, 13 high-order compact difference techniques have been included in more recent DNS 
efforts. Spectral-element methods have been used for spatial discretization around com- 
plex geometries. 14,15 Chebyshev collocation techniques have been used in boundary-layer 
and channel flow problems, 16 and a fully 3D non-periodic code has used used Chebyshev 
collocation techniques in two coordinate directions and compact and finite differences 
in the third coordinate direction. 17 Numerous DNS studies have used Adams-Bashforth, 
Runge-Kutta, Crank-Nicolson, and time-splitting approaches for time advancement. 

In the future, more advanced approaches should be considered. For example, a multi- 
grid general geometry DNS approach 18-21 has been developed and tested on a workstation. 
This multigrid DNS (MGDNS) approach has enabled computations of transition induced 
by 2D and 3D roughness and 3D transition on an unswept airfoil. 

Benefits obtained through careful programming and by updating the algorithms of a 
DNS code can be demonstrated by the example in Table 3. Using a Cray YMP super- 
computer, a DNS code was updated over a five-year period. Using the same grid and test 
problem, Table 3 clearly shows that both memory reductions and speed enhancements 
can be realized with minimal effort. (Note, the increase in cost per time step is attributed 
to a more expensive time-step procedure which buys a significantly lower number of time 
steps per period of computation.) 





Table 3. Program and Algorithmic developments for use in DNS. 



1990 

1995 

Memory (Mb) 

87.6 

35.7 

Startup Cost (sec) 

1896 

51 

Cost/Time step (sec) 

20.3 

48.37 

Cost/Period (hr) 

17 
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Often Poisson or Helmholtz equations must be solved during the course of a DNS. A 
Gauss-Seidel-like line iteration procedure 22 and direct solvers 17,23 have been used in DNS 
codes for Poisson/Helmholtz equations. With engineering and scientific applications as the 
customer, fast serial and parallel high-order direct solvers for Poisson 24 and Helmholtz 25 
equations (Dirichlet and Newmann boundary conditions) have been tested for speed and 
accuracy. 

Using computer-specific libraries and attributes can buy DNS significant improve- 
ments in compute power. For example, the CRAY BUFFER IN / OUT allows asynchronous 
read/ writes to disk. When the command is executed, data on a disk can be accessed si- 
multaneously while the code “crunches” numbers. On parallel computers, asynchronous 
communication enables computations to proceed while information is exchanged among 
the nodes. 

6. BOUNDARY CONDITIONS 

Boundary conditions have always been a critical issue in the use of DNS. Kloker et 
al. 22 summarized most of the available outflow treatments: outflow moving boundaries 26 
and the imposition of a non-physical region at the end of the physical region. In 
this non-physical region (or zone), a favorable pressure gradient or suction, 26 increased 
viscosity, 19 ’ 21 the Fringe method, 27 and the buffer domain technique 28 have been demon- 
strated to minimize non-physical reflections at the outflow. 

Pruett et al. 5 noted the boundary condition difficulties with high speed applications 
due to artificial Mach waves (and reflections) if the disturbance forcing, far-field boundary 
conditions and downstream conditions are not carefully handled. 

In the far-field, the disturbances are generally assumed to vanish so either homoge- 
neous Dirichlet or exponentially decaying boundary conditions have been used. But these 
assumptions can lead to considerable errors when the nonlinear effects are large and the 
mean-flow distortion quantity is important. 29 

7 FUNDAMENTAL FLUID MECHANICS 

There are many issues related to fundamentals in fluid mechanics where DNS (and 
LES) can play an important role. For example, we know that the flow along the 
attachment-line of swept wings is laminar or turbulent depending on the Reynolds num- 
ber. The experimental and theoretical results agree for the linear critical Reynolds 
number where small-amplitude disturbances become unstable on the attachment line at 
Rg — 245. 30-34 For large-amplitude disturbances, turbulence decays below some criti- 
cal Reynolds number and transition to turbulence will occur above this point. At this 
critical point, transition bypasses the conventional linear instability breakdown process. 
The experiments 35-38 showed that disturbances are damped for Rg < 100 and the flow 
becomes turbulent for Rg > 100. Note, the wide gap between the linear critical Reynolds 
number of Rg ~ 245 and the turbulent suppression critical Reynolds number of Rg ~ 100. 
Bridging this gap is important for wing design. Because of the wide spectrum of potential 
disturbances in the flow, DNS is required to resolve this gap. 




Concerning the fluid mechanics phenomena of local/global and absolute/convective 
instability, 39 DNS can play an important role in validating new theoretical concepts where 
traditional temporal or spatial theoretical approaches are alone invalid. For example, 
simulations which excite (forced actuation) global instabilities (in regions where the global 
instabilities may be marginally stable) can aid enhanced mixing. 

The unsteady flow about complex configurations have been studied 14,15 using DNS 
and should continue toward understanding physics and to devise/ validate simplified design 
tools. For example, impulsive flows 14 (transient phenomena) have been computed about 
a complex wedge configuration. 

7.1 ACTIVE FLOW CONTROL 

In the last decade, increased attention has been devoted to the development of me- 
chanical, distributed, and micro-sized techniques capable of enhancing our ability to mea- 
sure and control the unsteady flow in a wide variety of configurations (e.g., engine inlets 
and nozzles, combustors, automobiles, aircraft, and marine vehicles). Controlling the flow 
in these configurations may lead to greatly improved efficiency and performance, while 
decreasing the noise levels generally associated with the otherwise unattended unsteady 
flow. Depending on the desired result, one might wish to delay or accelerate transition, re- 
duce drag or enhance mixing. Furthermore, high-performance aircraft maneuvers without 
conventional mechanical devices may be a future goal (DOD, NASA, and industry funded 
programs). Removing conventional controls would lead to significant weight reductions. 


Blowing and/or suction 



Figure 3. Sketch of fluidic and surface displacement actuators. 

7.1. Unsteady Actuator-Induced Flow Control 

Aerodynamic design has a whole new set of challenges with the introduction of un- 
steady flow control. In particular, DNS (and LES) are required to identify the flow physics 
governing the coupled actuator-induced and shear flows (laminar through turbulent flow). 
The unsteady fluidic and surface displacement actuators (see figure 3) may involve zero, 
positive, or negative net mass, range in size from MEMS (less than 1 centimeter) to macro 
devices (centimeters), and operate at frequencies from 0 Hz to mili-Hz to kHz. Such active 
flow control concepts have been experimentally tested on an airfoil in a wind-t unn el . 40,41 
The results suggest that enhanced lift over the angle-of-attack range is available through 
the use of unsteady fluidic actuators. 

DNS can certainly be used to guide such experiments and to explain the complex 
fluid dynamics of the actuated flow field. To demonstrate the value of DNS, unsteady 



incompressible Navier- Stokes simulations were accomplished for the flow field induced by 
a single zero-net-mass actuator. 42 The actuator had a distribution of 

u(t,y )\ x = o = U m sm(ut) -sm(y/d) (1) 

where U m is the maximum actuator velocity, ui is the forcing frequency, and d is the 
actuator orifice size. Here, the forcing distribution is simply half of a sine-wave period 
in the y-direction; however, the simulations can be repeated when experimental data 
suggests an an alternate distribution. For this example problem, U m = 20m/ s, u = 0.5 
(~ 3183 Hz), and d = 0.5 mm are selected since these parameters are close to the family 
of synthetic-jet actuators. 43 

The computational domain had 201 grid points covering 18 mm in the jet-flow di- 
rection and 81 collocation points covering 16 mm in the cross-jet direction. The buffer 
domain begins 13.5 mm away from the actuator. Figure 4 shows contours of the u and 
v velocity components after 12 periods of forcing. Consistent with the experiments, 43 a 
jet-like flow develops in the simulations in spite of the zero net mass entering the com- 
putational domain. After reaching a peak amplitude close to the orifice, the amplitude 
of the jet begins to decay. Interesting to this fluidic actuator, a fairly strong transient 
phenomena is observed (in the middle of the computational domain) during the start-up 
of the actuation. 



7.2 Control Theory 

The coupling of control theory with DNS (or other unsteady nonlinear CFD approach) 
is necessary to expand our understanding of what kind of actuation is required for a given 
objective function. Being that actuation must be positioned on a surface of the configura- 
tion, the computations with control would be activated with a boundary condition. Such 
a concept would take the following form over a desired time interval t € (T 0 ,Ti), where 
0 < To < Ti < T. 

, = <Si m(To.T,) 

“' lr * \0 in(0,T„) and(Ti,T) ' ’ 

Here, the control function gi(t,x ) gives the rate at which fluid is injected, sucked, or the 
surface is displacement through the actuator T a (prescribed dimensions of the actuator). 




The actuator response based on the flow and objective function must be understood as a 
product of the optimization process. 

For example, a self-contained methodology for active flow control 44,45 was developed 
using control theory and DNS. The method coupled the time-dependent Navier-Stokes 
system with the adjoint Navier-Stokes system and optimality conditions from which op- 
timal states (i.e., unsteady flow fields and controls) were determined. The problem of 
transition delay was used as the initial validation case to test the methodology. The 
objective of control was to match the wall-normal stress to the steady base flow at a 
prescribed sensor location on the wall. Unlike feedback control methodologies wherein 
the sensed data determines the control through a specified feedback law or controller; the 
time-dependence of the control was the natural result of the minimization of the objective 
functional. Penalty terms were added to the functional to limit the size of the control 
gi and the oscillations in time. Constants (for example, a{ and /?,) were used to adjust 
the relative importance of the terms appearing in the functional. A sensor was used to 
feed information to a controller that in turn fed information to the actuator. (In the 
optimal control setting, the sensor is actually an objective functional and the controller 
is a coupled system of partial differential equations that determine the control that does 
the best job of minimizing the objective functional.) Control was effected through the 
injection and suction of fluid through a single orifice on the boundary. 

The method of Lagrange multipliers was formally used to enforce the constraints. 
Each argument of the Lagrangian functional was considered to be an independent vari- 
able so that each may be varied independently. For the incompressible Navier-Stokes 
equations, u,- were Lagrange multipliers that were used to enforce the a;, -components of 
the momentum equations, p was a Lagrange multiplier that was used to enforce the con- 
tinuity equation, and Si are Lagrange multipliers that were used to enforce the individual 
components of the boundary conditions (eqn. 2). The first variations of the Lagrangian 
with respect to the state variables u, u, and p were set equal to zero, leading to the adjoint 
or co-state equations. 


, du. 


dt 3 dxi 3 dxj dxi dxjdxj ’ dx t 
where the equations hold for t € (0,T). 


dui dp 


Ui m 
*lt=T 


were initial conditions and homogeneous conditions were used everywhere on the compu- 
tational boundary except where the sensor responds to the adjoint system. In the sample 
problem, the following forcing condition was used for the adjoint system. 

ui = ai ( n - Ti ) on ( T a ,T b ) x , (5) 

Ti(t,x ) is a given function (desired result) defined on ( T a ,Tb ) x T s , where (T a ,Tj) is a 
time interval such that 0 < T a < Tb <T. 

The optimality conditions were determined in a manner similar to the adjoint equation 
and are described by the equations below. 


— faT + 9i = ~J i y i ) on ( r o,Ti)xr a 

and gi(t,x ) is subject to initial and boundary conditions 

„ i _ _ _ /„\ _ i dgi\ 


9i\t=T 0 = 9i,o(x ) and 


dt \t=T 1 


(7) 



Given t at each point x 6 T a , a two-point boundary value problem in time over the 
interval (To, Ti). 

The coupled system was solved in an iterative manner. First, the simulation started 
with no control (<7; = 0) and the Navier-Stokes equations were solved for the velocity (u, v ) 
and pressure ( p ) fields. The adjoint equations were then solved for the co-state variables 
(u,v) and p. Then, using these adjoint variables, the optimality conditions were solved. 
The procedure was repeated until satisfactory convergence is achieved. 

Here, the optimal control methodology was applied to a boundary-layer flow having 
a single instability wave that can be characterized by a discrete frequency within the 
spectrum. These discrete small-amplitude instabilities can be suppressed through wave 
cancellation (WC), or wave superposition, using known exact information concerning the 
wave. Hence, the optimal control is “known” and can be used to validate the present 
DNS/optimal control theory approach in which the instability is to be suppressed with- 
out any a priori knowledge of said instability. For the computations, the grid has 401 
streamwise and 41 wall-normal points. The free-stream boundary is located 75 8* from 
the wall, and the streamwise length is 224<5* which is equal to approximately 8 Tollmien- 
Schlichting (TS) wavelengths. The nondimensional frequency for the forced disturbance 
is F = ui/R x 10 6 = 86; the forcing amplitude is vj = 0.1%. The Reynolds number based 
on the inflow displacement thickness (6*) is R = 900. To complete one iteration of the 
active-control simulation process, 4.5 min. on a Cray C-90 are required using a single 
processor for a iteration defined as three periods (T a — ► T& = 3T p ) of disturbance forcing. 

Based on 012 = 1 (selected arbitrarily) and #2 = 16.5 (refined through feedback 
method), the converged solution was obtained in 4 iterations. The wall-normal velocity 
(v) at a fixed distance from the wall and the actuator response (<72) is shown in figure 5 and 
compared with the desired wave-cancellation (WC) solution and a amplified uncontrolled 
solution. These results yield a delay in the transition by-way-of a suppression of the 
instability evolution. 




Figure 5. Disturbance wall-normal velocity and actuator response with downstream dis- 
tance for control (C3) and wave cancellation (WC). 44 (F-forcing location, A-actuator 
location, and S-sensor location). 





The ultimate goal of this line of research is to introduce automated control to ex- 
ternal flows over realistic configurations such as wings and fuselages, and to guide the 
development of actuators and control laws. The methodology may potentially be applied 
to separation control, re-laminarization, and turbulence control applications using one 
or more sensors and actuators. As DNS becomes useful in larger, more complex prob- 
lems (by-way-of computer speed advancements, parallel computing, advanced algorithm 
development), such controller approaches can guide active flow control projects. 

Finally, reducing the memory requirements of such methodologies is a primary task 
for future research. The adjoint system (eqn. 3) requires that the velocity field (it, v) 
obtained from the Navier-Stokes equations be known in the computational domain for all 
time. For this sample 2D problem, the iteration sequence and a modestly course grid, 246 
Mb of disk (or runtime) space are required to store the velocities at all time steps and for 
all grid points. Clearly for 3D problems the memory requirements becomes prohibitively 
expensive. This limitation can easily be removed if the flow-control problem involves 
small- amplitude unsteadiness (or instabilities). The time-dependent coefficients of the 
adjoint system reduce to the steady-state solution and no addition memory is required 
over the Navier-Stokes system. Perhaps, the coefficients may be stored every 10 (or more) 
time-steps or may be represented by a statistical average of the flow, thereby reducing the 
memory requirements by an order of magnitude. This hypothesis will require validation 
in a future study. 

CONCLUSIONS 

Based on the author’s subjective view, DNS will have limited growth on one-processor 
computers; however, parallel computers will enable solutions to a larger class of problems. 
Nonconventional methodologies (e.g., multi-grid) will lead to DNS solutions at an afford- 
able cost. DNS must be used in problems where simplification to the governing unsteady, 
nonlinear equations have not as yet been adequately validated. In the 1990s and beyond, 
such flows include active flow control on simple and complex configurations and traditional 
transition and turbulent shear flows on complex configurations. 
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